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FIB  STATUS  REPORT 


Tlie  attached  report*  constitutes  the  most  detailed  description  that 
we  have  produced  to  date  of  our 

Fields  by  Information  Blending  (FIB)  Methodology. 

We  believe  this  new,  comprehensive  approach  to  be  the  most  advanced 
and  definitive  available  for  the  objective  analysis  of  distributions  of 
specific  types  of  object  parameters.  The  methodology  is  applicable  to 
any  type  of  object  parameter,  scalar  or  vector,  having  significant  character 
expressible  by  linear  spatial  operators  acting  on  the  object  distribution 
(e.g. ,  information  as  to  gradient  components,  Laplaclan  and  other  higher 
derivatives,  divergence,  curl,  etc.). 

The  methodology  is  based  on  the  treatment  of  Information  as  a 
metered  commodity:  Every  input  piece  of  information  is  weighted  as  to  its 
purported  independent  worth.  This  reliability  weight  is  defined  as  the 
inverse  of  the  variance  associated  with  the  piece  of  information.  Information 
in  direct  observations  (i.e. ,  spot  measurements)  of  the  object  parameter, 
and  information,  observed  and/or  constructed,  as  to  spatial  characteristics 
of  the  object  distribution,  are  assembled  by  reliability  weights  and  blended 
to  give  the  best  resultant  distribution  of  the  object  value,  and  the  associated 
distribution  of  the  remaining,  unresolved  variance. 

The  attached  report  deals  with  the  application  of  the  methodology 
in  the  development  of  an  objective  analysis  scheme  for  Sea-Level -Pressure 
(SLP)  distribution.  References  are  made  to  the  simpler  application  of  the 


The  attached  report  has  also  been  reproduced  under  the  covers  of  the  Fleet 
Numerical  Weather  Central,  as  Technical  Note  No.  72-2. 


methodology  for  Sea-Surface-Temperature  (SST)  analysis*,  and  to  the  more 
complex  application  for  horizontal  wind  (UV)  analysis. 


At  this  time,  the  following  applications  are  in  operational  use  by 
the  Fleet  Numerical  Weather  Central: 

FIB/SST  Northern  Hemisphere  Oceans,  125x125  grid  array  (half  mesh 
length) 

FIB/SST  Gulf  Stream  Region,  63x63  grid  array  (quarter  mesh  length) 

FIB/SLP  Northern  Hemisphere,  125x125  grid  array  (half  mesh  length) 

The  mesh  length  refers  to  the  grid  spacing  of  the  standard  FNWC  63x63  grid 
array  of  the  polar  stereographic  projection  encompassing  the  northern  or 
southern  inscribed  hemisphere. 

The  following  applications  are  currently  under  development: 

FIB/SST  Southern  Hemisphere  Oceans,  63x63  grid  array  (standard  mesh 
length) 

FIB/OTS  Analysis  of  Ocean  Thermal  Structure  (OTS)  Parameter  Distributions 

FIB/SLP  Med.  Greater  Mediterranean  Region,  63x63  grid  array  (quarter 
mesh  length) 

following  applications  have  been  formulated: 

Horizontal  Wind  Analysis 

Upper  Air  Analysis,  isobaric  height  and  horizontal  wind 
structure  parameters 

Other  applications  are  under  formulation . 


The 

FIB/UV 

FIB/UA 


Holl,  Manfred  M.,  Bruce  R.  Mendenhall  and  Charles  E.  Tilden,  1971; 
"Technical  Developments  for  Operational  Sea  Surface  Temperature  Analysis 
with  Capability  for  Satellite  Data  Input"  ,  Final  Report,  Contract  No. 
N62306-70-C-0334 ,  (Naval  Weapons  Engineering  Support  Activity), 
Meteorology  International  Incorporated,  Monterey,  California. 
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1. 


Introduction 


A  comprehensive  technique  for  the  objective  analysis  of  scalar  and  of 
vector  fields  has  been  developed  by  Meteorology  International  Incorporated, 
under  Navy  contracts.  This  technique  has  been  designated  the  Fields  by 
Information  Blending  (FIB)  technique . 

Information  is  a  delicate  commodity  which  is  easily  depreciated  by 
transformation  and  assimilation  schemes.  The  FIB  technique  treats  information 
as  a  metered  commodity.  It  assembles  all  available  information,  which  must 
include  suitable  weighting,  and  blends  the  information  into  resultant  analyses 
and  associated  distributions  of  resultant  resolution  weight. 

In  the  FIB  context,  all  information  statements  must  include  both 
parameter  estimates  and  associated  reliabilities.  Any  statement  of 
information,  whether  observation  or  analysis  value  at  a  grid  point,  is 
incomplete  without  an  associated  reliability. 

For  an  independent  piece  of  information  the  reliability,  or  report 
weight,  is  defined  as  the  inverse  of  the  error  variance  inherent  in  the 
observation  and/or  associated  with  the  class  of  observation.  For  a 
statement  of  resolution  the  associated  reliability,  or  resolution  weight, 
is  the  inverse  of  the  unresolved  variance. 

The  FIB  technique  is  based  entirely  on  rules  for  adding  uncorrelated 
variance  contributions  and  for  adding  independent  information.  The 
operation  which  produces  the  implied  resultant,  given  weighted  grid-point 
arrays  of  values  of  the  object  parameter  and  of  its  linear  integral  and 
finite-difference  derivatives,  is  called  the  blending. 

The  information  available  for  producing  an  analysis  may  be 
considered  to  be  from  three  sources  in  time: 

(1)  Concurrent  observations.  Observations  taken  at  or  very  near 
the  applicable  time  of  the  object  analysis. 
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(2)  Near- past  observations.  This  information  has  been  assimilated 
into  the  analysis  preceding,  and  is  projected  forward  along  the 
time  axis  in  the  form  of  first-guess  fields. 

(3)  Well-past  observations.  This  is  the  source  of  physical  and 
statistical  relationships .  A  pertinent  example  is  the  geostrophic 
approximation  by  which  wind  information  is  transformed  into 
pressure  gradient  information  and  vice  versa.  Note  that  the 
transformation  contributes  an  error  variance  to  the  information 
so  transformed.  Well-past  observations  are  also  the  source  of 
specific  inferences  such  as  those  grouped  under  area  synoptic 
skills  of  an  analyst. 

The  information  available  also  differs  as  to  type  and  relevance.  The 
FIB  technique  can  directly  assimilate  information  in  observations  (i.e. , 
random  samplings)  of  the  parameter  and  information  as  to  integral  and 
differential  field  properties. 

Two  versions  of  the  FIB  technique  have  been  accomplished  to  date: 

FIB/SST:  Application  to  Sea-Surface  Temperature  with  a  hemispheric 

(125x125  grid  array)  and  a  fine  mesh  (63x63)  regional  adaptation. 

FIB/SLP:  Application  to  Sea-Level  Pressure  with  a  hemispheric  (125x125) 
adaptation . 

A  third  version  has  been  formulated: 

FIB/UV:  Application  to  a  horizontal-velocity  field  (i.e. ,  wind  analysis) . 

FIB/SST  directly  assimilates  information  in  the  form  of  (1)  temperature 
values,  and  (2)  temperature-gradient  values  in  two  components.  The 
temperature  gradient  is  expressed  in  terms  of  temperature  difference  per 
grid  interval.  It  has  two  components,  associated  with  the  two  grid-array 
coordinate  axes. 


-2- 


FIB/SLP  directly  assimilates  information  in  the  form  of  (1)  pressure 
values,  (2)  pressure-gradient  values  in  two  components,  such  as  derived 
from  winds,  and  (3)  Laplacian  values  for  pressure,  expressed  in  terms  of  a 
five-point  finite-difference  operator,  in  the  grid  scale. 

FIB/UV  directly  assimilates  information  in  the  form  of  (1)  wind  values 
in  two  components,  (2)  vorticity  in  the  grid  scale,  and  (3)  divergence  in  the 
grid  scale. 

Figures  1,  2  and  3  relate  to  FIB/SST,  FIB/SLP  and  FIB/UV  respectively. 
Each  indicates  the  reference  locations  for  parameter  values  in  an  arbitrary 
area  module  of  a  rectangular  grid  array. 

The  FIB  technique  includes  the  capability  to  accept  artificial  data 
whether  generated  objectively  or  subjectively  by  interpretation  of  additional 
information  such  as  from  satellites  or  by  preferred  interpretation  of  surface 
data.  The  artificial  data  must  be  expressed  in  the  direct  forms  and  must 
include  estimated  reliability  weights . 

In  the  blending  operation  the  spreading  of  observations  (i.e. , 
individual  independent  reports)  is  effected  by  the  highest-order  fields — that 
is,  by  the  resolution  fields  which  express  information  derived  from  higher- 
order  differential  properties  of  the  field  which  control  shape  and  character: 

For  FIB/SST  this  is  the  gradient  resolution  in  two  components. 

For  FIB/SLP  this  is  the  Laplacian  resolution  field. 

For  FIB/UV  this  is  expressed  by  the  vorticity  resolution  field  and  the 
divergence  resolution  field. 

These  are  the  basic  spreaders.  The  blending  of  observations,  and  the 
exploitation  of  first-guess  information,  can  be  enhanced  by  expressing 
information  in  additional  differential  properties  of  the  fields. 
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The  FIB  technique,  in  general,  includes  the  following  component 
operations: 

(1)  Preparation  of  first-guess  fields  from  available  and  generated 
sources . 

(2)  Assembly  of  new  information.  Independent  estimates  in  any 

of  the  acceptable  forms  are  assembled  by  weighted  combination 
at  grid-point  reference  locations,  but  the  highest-order  fields 
must  be  restructured  to  accept  new  information. 

(3)  Blending  of  information  fields  into  resultant  analysis. 

(4)  Evaluation  of  reliability  weights,  performed  by  limited 
reanalysis  using  a  perturbed  value  at  each  grid  point. 

(5)  Gross-error  checking  and  reevaluations.  Includes  capability 
for  monitoring  quality  of  station  reports . 

(6)  Recycling  of  Assembly,  Blending  and  Resultant-Weight 
operations. 

These  operations  are  detailed  for  FIB/SST  in  a  technical  report*.  The 
operations  for  FIB/SLP  are  detailed  in  the  present  report. 

A  major  source  of  information  for  sea-level  pressure  analysis  is 
represented  by  surface  wind  observations.  FIB/SLP  exploits  much  of  this 
information  by  transforming  each  wind  into  the  two  pressure-gradient 
components,  according  to  balanced  large-scale  dynamics,  neglecting  local 
wind  tendency.  Those  surface  winds  which  appear  to  be  non-representative 
of  the  large-scale  dynamics  may  be  rejected  in  the  gross-error  check 
operation.  In  the  synoptic  scale  this  treatment  exploits  most  of  the 
relevant  information  inherent  in  the  winds. 


*Holl,  Manfred  M.  ,  Bruce  R.  Mendenhall  and  Charles  E.  Tilden,  1971; 
"Technical  developments  for  operational  sea  surface  temperature  analysis 
with  capability  for  satellite  data  input".  Final  Report,  Contract  No. 
N62306-70-C-0334  (Naval  Weapons  Engineering  Support  Activity) , 
Meteorology  International  Incorporated,  Monterey,  California,  46pp. 
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FIB/SST 


Fig.  1  The  arbitrary  Z,m  area  module  (consisting  of  one  comer  point, 

two  sides,  and  the  interior  area)  and  adjoining  modules.  Reference 
locations  in  the  module  are  shown  for  the  Sea-Surface-Temperature 
parameter,  T^  ^  ,  and  the  finite-difference  parameters: 


T  -  T 

l,m+\  Z,m 


and 


T 

SL+l  ,m 


The  respective  reliabilities  of  T,  |i  and  v  values  are  denoted  by 
A,  B  and  C;  the  reliability  of  a  value  is  defined  as  the  inverse  of 
the  error  variance  associated  with  the  value. 
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FIB/SLP 


Fig.  2  The  arbitrary  i  ,m  area  module  (consisting  of  one  comer  point, 

two  sides,  and  the  interior  area)  and  adjoining  modules.  Reference 
locations  in  the  module  are  shown  for  the  Sea-Level-Pressure 
parameter,  P.  ,  and  the  finite-difference  parameters: 


l  ,m 

=>  p 

j l  ,m  +  \ 

-  P  '  v  => 

i ,m  '  i  ,m 

P  -  P 

i+\  ,m  SL ,m  ' 

and 

L . 
l  ,m 

P  +  P 

JL+l,m  &  ,m+ 1 

+  P .  .  +P.  ,  ■  4  P. 

i,m- 1  i,m 

The  respective  reliabilities  of  P  ,  \x  ,  v ,  and  L  values  are  denoted 
by  A,  B,  C  and  D;  the  reliability  of  a  value  is  defined  as  the 
inverse  of  the  error  variance  associated  with  the  value. 


-6- 


FIB/UV 


Fig.  3  The  arbitrary  Z,m  area  module  (consisting  of  one  comer  point,  two 
sides,  and  the  interior  area)  and  adjoining  modules.  Reference 
locations  in  the  module  are  shown  for  the  u  and  v  components  of 
the  wind  velocity,  and  the  finite-difference  parameters: 

e  =>  u  —  u  +v  —  v  : 

Z  ,m  Z+l,m  l  ,m  Z,m  + 1  Z  ,m 


and 


Z,m 


v 


Z,m 


v 


Z-\,m 


u 


Z,m 


u 


Z  ,m  -1 


The  respective  reliabilities  of  u,  v,  q  and  e  values  are  denoted 
by  B,  C,  D  and  E;  the  reliability  of  a  value  is  defined  as  the 
inverse  of  the  error  variance  associated  with  the  value. 
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2. 

Conceptions  and  Basic  Formulas 

2.1 

Definitions  and  Rules 

2.1.1 

Definitions 

Let  p  denote  an  independent  estimate  (e.g.  ,  measurement)  of  sea- 
n 

level  pressure.  Let  o  denote  the  standard  deviation  associated  with  this 

n  2 

estimate.  The  variance  is  given  by  .  The  weight  or  reliability  of  the 
information  expressed  by  the  estimate  is  defined  as  the  inverse  of  the 
variance: 


A  =  1/ct  2 
n  n 


(1) 


Let  u  denote  an  estimate  of  the  finite  difference  in  pressure  between 
m 

two  locations,  m  and  m  +  1 : 

u  =»  an  estimate  of  p  ,  ,  -  p  .  (2) 

m  m  + 1  m 


The  weight  of  a  is  denoted  by  B  .  The  associated  variance  is  then  given 

_  171  171 

by  B  '  . 

m 


2.1.2  Addition  of  Contributing  Variances 


Estimates  are  independent  if  the  associated  errors  are  mutually 
uncorrelated.  For  a  linear  combination  of  such  independent  estimates, 

P  =  PQ  ±  M !  ±  M  2  —  •  •  •  /  (3) 
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the  associated  variance  is  given  by 


Or2  =  A_1  =  Ao_1  +  B^1  +  B2_1  +  .  .  .  .  (4) 

For  the  product  of  an  estimate,  p  ,  multiplied  by  a  constant,  c  ,  the 

2  2  ’■•••'  ••  ‘ 

variance  is  c  cr  . 


2.1.3  Relevant  Variance  in  a  Measurement 


With  respect  to  the  analysis  of  a  spatial  distribution,  the  contributing 
variances  in  a  measure,  relative  to  an  objective  resolution  as  to  scale,  are  . 

-  V- "  •  •  >'  ■  ■  ••• 

2  2  2 

a  -  a  (instrument  system)  +  a  (subscale  noise)  .  (5) 

Should  the  scale  of  the  objective  resolution  be  extended  in  detail 
(i.e. ,  in  smallness) — as  made  possible,  for  example,  by  a  finer  grid  array — 
then  the  noise  component  reduces  and  the  relevant  information  weight  of  the 
measure  is  increased. 


2.1.4  Addition  of  Information 

Let  p  ,p2,.  .  .p^,.  •  •Pjj  be  N  independent  estimates  of  an 
object  value  p  .  Let  the  respective  weights  be  A^A^.  .  .A^, .  .  .A^  . 
The  best  estimate  of  p  that  this  information  affords  is  expressed  by 


A^i  +  ^2^2 


A  p  +  .  . 
n  n 


ANPN 


ai  +  a2  + 


.  A 


N 


.  .  A  + 
n 


(6) 


And  the  weight  is  expressed  by 


A  =  A  +  A  +  .  .  .  A  +  .  .  .  A 
.1.2  n  N 


This  assimilation  can  be  performed  sequentially,  saving  only  the  accrued 
value  and  its  growing  weight.  The  result  is  independent  of  the  order  of 
assimilation. 


2.1.5  Removal  of  Information 


The  contribution  of  a  piece  of  information ,  p  of  weight  A  , 

n  n 

which  has  already  been  assembled  into  a  resultant,  p*  of  weight  A*,  can 


be  removed  by  subtraction: 


A *p*  -  A  p 
_ n  n 

A*  -  A 


A*  -  A 


The  information,  p  of  weight  A  ,  represents  the  independent  background 

B  B 

(i.e. ,  verifying)  information  relative  to  the  nth  report  which  had  been 
assembled  into  the  resultant. 

For  a  report  which  was  withheld  from  the  assembly  the  total  resultant, 
p*  of  weight  A*  ,  represents  the  independent  background  information: 


PB  =  P 


A  =  A 

ab  a  ' 


The  background  information  may  be  used  to  gross-error  check  and 
reevaluate  any  piece  of  information. 
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2.2 


Gross-Error  Check  and  Reevaluation 


A  piece  of  information  is  invalid  if  it  is  not  what  it  purports  to  be. 

If  the  reported  value  differs  widely  from  that  measured,  or  if  the  reported 
location  of  the  measurement  is  considerably  removed  from  the  true  one, 
then  the  report  is  termed  a  gross  error.  For  an  objective  analysis  technique 
to  be  acceptable  it  must  have  a  well-developed  capability  for  invalidating 
gross  errors. 

A  piece  of  information  which  lies  in  the  wings  of  the  distribution 
implied  by  the  associated  variance  is  termed  an  outlier  or  maverick.  It  is 
a  valid  report  but  is  not  representative. 

Gross  errors  may  occur  with  any  value.  They  are  distinguishable 
only  if  the  disparity  with  the  background  information  is  significant.  Mavericks 
are  always  confused  to  some  extent  with  gross  errors.. 

Let  p  of  weight  A  be  a  piece  of  information  (e.g.  ,  the  nth  report) 
n  n 

which  is  to  be  checked  for  validity  and  reevaluation  (confirmation  or  weight 
reduction).  Let  p  of  weight  A  be  the  independent  background  information 

D  D 

relative  to  the  same  object — but  unknown — value  p  .  The  actual  difference 
between  the  estimates  is 

pn  '  PB  •  .  (11) 

2 

The  associated  variance,  here  denoted  by  cr  ,  of  this  difference  is 

n ,  d 


n 


+  A 


-1 

B 


A_  +  A 
B  n 

A  A 
n  B 


(12) 
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The  actual  difference  measured  in  units  of  the  associated  standard  difference, 
£T  _  ,  is  denoted  by  X  .  This  normalized  difference  is  given  by 


X 


2 

n 


A 

n 


(13) 


The  purported  weight,  A  ,  is  reevaluated  on  the  basis  of  X 

n  n 

The  reevaluated  weight  is  expressed 


A  D  =  ( 

nR 


1  +  X 


A  for  X  =  1 
n  n 


-  A  for  1<X2^X2 
2  n  n 

n 


2  2 

0  for  X  >  X 
n 


(14) 


2 

For  X  ^  1  the  nth  report  is  deemed  to  be  acceptable  as  purported, 
n 

2  ~2 

For  1  <  X  <  X  the  nth  report  is  deemed  to  be  valid,  but  in  order  to 
n 

reduce  the  influence  of  mavericks  the  weight  is  reduced. 

2  ~2 

For  X  >  X  the  nth  report  is  deemed  to  be  a  gross  error  and  is  rejected 
n 

as  worthless. 

~~2 

The  specification  of  X  is  somewhat  arbitrary.  The  specification 

2 

is  based,  in  application,  on  examinations  of  distributions  of  X  for  sample 
populations  of  reports.  This  procedure  is  referred  to  as  tuning  the 
parameter. 
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2.3 


Propagation  of  Information 


2.3.1  Explicit  Blending 


Values  Weights 


1  o 

1 

1 

P1 

Ai 

1 

1 

1 

H  * 

P2  "  P1 

Bi 

1 

2  o 

i 

i 

P2 

A2 

i 

i 

i 

^2  " 

P3'P2 

B2  , 

i 

3  o 

i 

P3 

A3 

i 

i 

i 

^3  " 

P4  "  P3 

,  B3 

i 

4  o 

P4 

A4  ' 

• 

Fig.  4  A  one-dimensional  array  of  independent  information 


Consider  the  implications  of  the  distribution  of  independent 
information  shown  in  Fig.  4.  The  pressure  information  at  any  one  point 
generally  has  inferences  at  the  other  points  via  information  in  the  form  of 
estimates,  \i  ,  of  the  pressure  difference  between  successive  points. 

The  estimates,  \i,  are  finite-difference  expressions  of  the  first  derivative 
of  pressure  along  the  axis  of  the  array. 

What  is  the  resultant  pressure,  p*  ,  and  the  resultant  weight.  A*  , 
say  at  level  4 ,  that  is  implied  by  the  available  information  ? 
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Four  estimates  are  available  for  p 


Value 


Weight 


Pl  +  "l  +  '*2  +  '13 


P2  +  P2  +  »I3 


P3  +  t*3 


Ar1+Bi1+B2'1+B3'1 


A2'1  +  B2'1  +  B3'1 


A3'1  +  B3_1 


However  it  would  be  wrong  to  combine  these  by  the  rule  for  adding  independent 
estimates — because  they  are  not  independent.  The  errors  are  correlated;  there 
are  common  contributions  to  the  variances  of  the  first  three  estimates  listed. 

There  is  a  stepwise  procedure  for  combining  the  available  information 
which  satisfies  the  condition  of  independence.  The  information  which 
propagates  toward  point  4  can  be  accrued  step  by  step.  At  point  2  the  values 


Pl  +  /x  i  of  weight  1  +  B.^  1 


p^  of  weight 


can  be  combined  by  the  rule  for  addition  of  independent  estimates: 


2(1+2) 


(£i£)  (pi+pi)  +  A2p: 

(Ar1+Bi_1) 1  +  a2 


2(1+2) 
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The  subscript  parentheses  have  been  added  to  enclose  the  sequence  of 
combination. 

Next,  at  point  3,  the  values 


p2(l+2)  +li2  o£wel9ht  (A 


-1 

2(1+2) 


+  B 


-1 


-1 


and 


p  of  weight  A 

<3  <3 


can  be  combined  to  form 


P3( 1+2+3)  ofwelght  A3  (1+2+3) 
The  procedure  is  repeated  to  form 


P4(l+2+3+4)  ofWeiSht  A4(l+2+3+4)  ' 

which  are  the  desired  resultants,  p*  and  A*  . 

The  information  which  arrives  at  point  4  from  the  one  direction  is 
the  value 


P3(l+2+3)  +  ^  3  °£  weight  (a'^+j,  +  B,"1) 


-1 


The  weight  is  upper- bounded  by  A3^+2+^  and,  more  significantly,  by  B3  . 
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For  an  intermediate  point,  say  point  2,  the  ambient  information 
arrives  from  two  directions.  The  available  independent  estimates  are 

P1  +  of  weight  1  +  Bj-1 

P3<4+3)-“2  weight  (a3(4+3)  +  B^1 ) 

P2  of  weight 

This  stepwise  procedure  can  be  applied  to  any  arrays  in  which  each 
point  is  connected  to  any  other  point  by  a  single  path.  The  resultant  at 
point  6  of  Fig.  5a,  for  example,  can  be  calculated  as  follows:  The 
information  from  points  1 ,  3+2 ,  and  4  are  combined  with  that  at  point  5 
for  an  accrual  value  and  weight  at  point  5 .  This  accrual  at  point  5  is  then 
combined,  at  point  6,  with  the  information  from  points  7  and  8,  and  that  at 
6,  giving  the  resultant  at  point  6. 

8 

(a)  (b) 

Fig.  5  Example  of  (a)  a  single-path  array,  and  (b)  a  multi- path  array. 

This  stepwise  procedure  is  termed  explicit  blending. 

The  explicit  blending  procedure  calculates  resultant  values  and 
resultant  weights. 

Explicit  blending  cannot  be  applied  to  multi-path  arrays . 
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2.3.2 


Related  Applications 


2.3.2. 1  Conceptual  Basis 

Explicit  blending  is  useful  for  combining  estimates  of  values  and 
first  differences  along  one  dimension.  It  has  found  application  in  merging 
new  data  with  climatology,  at  depths,  for  rendering  resultant  ocean  thermal 
structure*.  While  not  equal  to  the  present  primary  task  of  objective  analysis 
for  two-dimensional  arrays,  including  higher-order  difference  information, 
explicit  blending  forms  part  of  the  conceptual  basis.  It  also  arises  in 
several  related  considerations  in  sea-level  pressure  analysis. 


2. 3. 2. 2  Collections  of  Pressure  Reports  with  Reduction  to  Sea  Level 

A  related  consideration  arises  in  the  collection  of  pressure  reports 
in  the  realm  of  a  single  collecting  point  (i.e.  ,  a  grid  reference  location  for 
pressure)  in  regions  where  terrain  elevation  is  significant.  The  direct,  but 
inconsistent,  scheme  which  is  generally  followed  is  to  reduce  each  station 
pressure,  separately,  into  a  sea-level  pressure  report,  by  some  formula 
for  the  fictitious  pressure  addition  associated  with  the  drop  in  elevation  to 
sea  level.  Each  such  reduced  sea-level  pressure  report  has  a  variance 
made  up  of  the  station-report  variance  and  a  variance  addition  for  the 
reduction  to  sea  level.  The  translation  of  each  of  these  reduced  sea-level 
pressure  reports  to  the  nearest  grid  reference  location  for  pressure 
(See  Fig.  2) — using  available  information  as  to  the  local  shape  of  the 


*Mendenhall,  Bruce  R. ,  1970;  "Design  of  a  structure  parameterization  method 
for  application  to  the  three-dimensional  analysis  of  ocean  temperature".  Final 
Report,  Contract  No.  N00228-69-C-0833  ,  (Fleet  Numerical  Weather  Central) , 
Meteorology  International  Incorporated,  Monterey,  California,  57  pp. 
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sea-level  pressure  field — adds  another  contribution  of  variance.  The 
addition  of  the  collection  of  such  reduced  reports  by  the  rule  for  information 
addition  is,  generally,  in  violation  of  the  condition  of  independence  of 
estimates.  As  more  and  more  reports  are  collected,  the  accrued  weight  can 
grow  without  bound  even  though  the  minimum  terrain  elevation  may  be 
thousands  of  feet  in  the  region.  This  wrong  method  is  schematically 
depicted  by  Fig.  6a. 


(a)  (b) 


Fig.  6  Collection  of  Station  Pressures  at  Sea  Level 

The  correct  method  is  schematically  depicted  by  Fig.  6b.  Each 
station  pressure  should  first  be  projected,  at  station  elevation,  to  the 
nearest  grid  reference  location  in  the  horizontal.  This  projection  involves 
a  variance  contribution.  All  projected  reports  then  pertain  to  elevations 
along  the  vertical  at  the  grid  reference  location.  The  resultant  at  sea  level 
can  then  be  calculated  by  stepwise-downward  explicit  blending  to  sea  level. 
The  minimum  variance  in  the  resultant  is  that  associated  with  the  reduction 
from  the  lowest  station  elevation  to  sea  level;  no  matter  how  many  reports 
enter,  the  total  weight  is  thereby  bounded. 


-18- 


We  are  currently  using  the  inconsistent  scheme,  depicted  by 
Fig.  6a,  because  station  pressures  are  not  generally  transmitted.  They 
are  reduced  to  sea  level  at  the  station  for  transmission.  Formulas  used 
for  reduction  to  sea  level  differ  from  country  to  country  and  even  within 
countries.  Recovery  of  station  pressure  from  the  reduced  sea-level  pressure 
and  station  elevation  is  not  generally  feasible. 

We  are  considering  a  refinement  which  will  bound  the  total  collected 
weight.  According  to  this  refinement  each  station  report  weight  is  reduced 
by  adding  variance  associated  with  reduction  down  to  the  elevation  of  the 
lowest  station  in  the  collection;  variance  is  also  added  in  association  with 
the  horizontal  projection  to  the  grid  reference  location.  The  reports  so 
weighted  are  combined.  The  total  weight  is  then  reduced  by  adding  a 
variance  associated  with  the  remaining  reduction  to  sea  level.  Adaptation 
of  this  refinement  also  involves  some  modification  of  the  gross-error 
checking  and  reevaluation  procedures. 


2. 3. 2. 3  Collections  of  Wind  Reports 

Similar  considerations  arise  in  the  explicit  use  of  wind  observations. 
Each  wind  may  be  transformed  into  a  pressure  gradient  estimate  (two 
components)  via  a  balance  approximation.  The  associated  variance  is  a 
combination  of  that  equivalent  to  the  variance  in  the  wind  estimate  plus 
the  contribution  inherent  in  the  balance  approximation. 

Consider  a  collection  of  wind  reports  in  the  realm  of  a  single 
collecting  point  for  a  gradient  component.  There  is  a  wrong  way  and  a 
right  way  for  combining  the  information  into  a  gradient  estimate. 
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The  wrong  way  is  to  combine  the  transformed  winds  as  estimates  of 
the  gradient  component,  with  weights  which  have  each  been  reduced  for  the 
imbalance  contribution  of  variance.  The  more  winds  there  are  in  the  collection, 
the  greater  the  apparent  resolution  of  the  gradient  component,  for  the  weight 
accrues  without  bound.  This  is  unrealistic  because  there  must  remain  the 
unknown  imbalance  between  wind  and  pressure.  The  fault  lies  in  the  fact 
that  the  gradient  estimates  are  not  independent:  they  share  the  common 
variance  contribution  inherent  in  the  balance  approximation. 

.  The  consistent  way  to  combine  the  gradient  information  in  a  collection 
of  winds  is  as  follows:  Each  wind  is  transformed  by  the  balance  approximation 
into  an  equivalent  pressure  gradient,  but  the  variance  is  still  regarded  as  only 
that  contributed  by  the  wind  report  variance  now  expressed  in  the  equivalent 
pressure  gradient  units.  The  collection  of  winds,  expressed  as  equivalent 
pressure  gradients,  are  combined  with  the  weights  which  reflect  only 
equivalent  wind  variance.  This  gives  the  resolution  of  the  wind  expressed 
as  an  equivalent  pressure  gradient  with  weight  in  the  same  units.  The 
weight  is  unbounded  as  an  expression  of  wind  resolution.  The  last  step  is 
inherent  in  assuming  this  equivalent  pressure  gradient  to  be  an  estimate  of 
the  pressure  gradient.  This  assumption  involves  the  addition  of  the  variance 
inherent  in  the  balance  approximation .  The  resultant  weight  of  the  gradient 
estimate  is  thereby  upper  bounded. 

Strictly  considered,  the  wind  information  so  collected  at  one  grid 
reference  location  is  not  fully  independent  of  that  collected  at  a  nearby 
reference  location.  The  variance  contributed  by  the  imbalance  is  spatially 
correlated.  This  is  one  of  the  inconsistencies  we  have  accepted,  and  for 
which  partial  allowances  are  made. 
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2. 3. 2. 4  Collections  of  Off-Time  Reports 


The  design  of  schemes  for  using  off-time  reports  is  similarly  loaded 
with  inconsistency  traps . 

Consider  a  collection  of  off-time  pressure  reports  in  the  realm  of  a 
single  grid  reference  location.  Each  report  provides  an  estimate,  at 
analysis  time,  with  weight  reduced  by  the  addition  of  variance  growth  in 
the  time  disparity.  The  collection,  however,  cannot  be  combined  with  these 
reduced  weights,  because  the  off-time  variance  growth  is  correlated  between 
reports.  The  consistent  way  is  to  explicitly  blend  the  reports  along  the  time 
axis  at  the  grid  reference  location,  as  was  proposed  for  the  vertical  axis  in 
collecting  pressure  reports  with  reduction  to  sea  level,  in  Section  2. 3. 2. 2. 

Such  an  explicitly  blended  collection  of  off-time  reports  at  one 
reference  location  is  not  fully  independent  of  such  a  collection  at  a  nearby 
reference  location.  The  variance  growth  factor  is  spatially  correlated.  For 
short  time  periods  this  correlation  can  be  ignored.  As  the  period  approaches 
that  of  the  analysis  cycle  the  problem  transforms  into  that  of  generating 
first-guess  fields  (i.e. ,  information  carried  along  the  time  axis)  from  prior 
analyses  and  predictions  therefrom. 
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2.4 


Implicit  Blending 


2.4.1  Formulations 

The  FIB  technique  includes  an  assembly  operation,  followed  by  a 
blending  operation.  In  the  assembly  operation  each  new  piece  of  independent 
information  is  combined  with  other  independent  information  (including  first- 
guess  information)  at  the  nearest  reference  location  for  the  parameter.  The 
basic  parameters  and  reference  locations  are  shown  in  Fig.  2  for  the  Sea-Level 
Pressure  version.  The  array  of  assembly  locations  serves  as  a  repository  for 
the  independent  information  in  this  assembly  operation. 

The  analysis  is  effected  by  the  blending  operation.  The  object  of 
the  blending  operation  is  to  determine  the  best  resultant  field,  p*  ,  with 
associated  resultant  weight  field,  A*  .  The  concept  may  be  expressed  by 

A*  d*  —Ad  +  A3  D^ 

i ,m  vJL,m  i ,m  p i ,m  i ,m  pJL,m 

i  Independent  ambient  information 

i 

Independent  information  assembled  at  the  reference  location 


+ 


l  ,m 


(17) 


These  resultants  are  derived  at  every  reference  location  for  pressure.  The 
independent  ambient  information  stems  from  all  reference  locations  for 
pressure  other  than  1  ,m  and  is  propagated  by  gradient  and  higher-order 
difference  information. 

The  pressure  information  at  reference  locations  is  connected  by 
multiple  paths  and  multiple  finite-difference  parameters .  Explicit  blending 
is  not  equal  to  the  task  of  effecting  the  consistent  blending.  An  implicit 
blending  scheme  is  required. 


-22- 


Early  efforts*  to  formulate  a  consistent  implicit  blending  scheme 
proved  inadequate.  Later,  we  devised  and  tested  a  formulation  which  has 
stood  up  to  all  available  tests  for  consistency.  This  made  possible 
development  of  the  first  FIB  version,  that  for  sea-surface-temperature 
analysis**. 

The  formulation  is  based  on  minimizing  the  pertinent  error  functional — 
the  total  weighted  disparity  with  all  available  elements  of  information. 

The  basic  error  functional  for  sea-level-pressure  analysis  is  defined 
as 


The  summation  extends  over  all  area  modules.  The  elements  of  information 
are  referred  to  the  locations  shown  in  Fig.  2.  The  boundary  modules  are  not 
distinguished  except  that  weights  of  all  those  information  elements  which 
couple  with  the  exterior  are  identically  zero.  Each  objective  resultant 
element  of  p*  is  that  which  minimizes  its  contribution  to  E. 

*Tilden,  Charles  E.,  James  R.  Clark  and  Manfred  M.  Holl,  1969;  "Techniques 
for  the  utilization  of  satellite  observations  in  operational  mapping  of  sea 
surface  temperature.  Final  Report,  Contract  No.  N62306-C-0288  (Naval 
Air  Systems  Command),  Meteorology  International  Incorporated,  Monterey, 
California,  33  pp. 

**Ibid,  p.  4 
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2.4.2  The  Blending  System  of  Linear  Equations 

The  arbitrary  X  ,m  component  of  the  blending  system  of  linear 
equations  is  given  by  setting 


9E 


=  0 


dp 


X  ,m 


(19) 


This  gives 


S,  p, 

X ,m  X ,m 


AJL,m  pJL,m 


BX,m  (pX,m- 1  ~  ^Z,m)  ~  BX,/rc- 1  (pX,w- 1  +  ^X,m-l) 

^X,ot(PX+1,ot  V H,m)  ^X-l ,m  (pX-l,m  +  ^X-l,^) 


-  4  D 


(pJ. 


,  +  p*  ,  +  p*  ,  +  p*!  ,  -  L  „ 

w  +  1  *X+l,/rc  X,/rc-l  l-l,m  l  ,m 


) 


/  *  _  *  *  *  \ 
■  DX,m-l  rPX#w-l  +  Lx,w-1  "  PX+l,m-l  ”  PX,77j-2  ~ 


DX-1,ot  (4  +  Lx-i, 


m 


-  P 


+  l  "  pX-1,ot-1  ^X-2 ,m 


-  P, 


) 


DJi,m  +  l \4  pX,w  +  l  +  LX,m  +  l  "  PJL,m+2  "  pX+l,w+l  "  PX 


-l,OT  +  l) 


DX+l,w  (4  pX+l,m  +  LX+l,/rc  "  pX+l,/rc  +  l 


*  *  \ 
"  pX+2 ,m  "  pX+l,^-i/ 


=  0 


(20) 
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where 


St  =  A.  +  B,  +  B.  ,  +  C,  +  C„  , 
&,m  H  ,m  i,m  JL,m- 1  JL-l  ,m 


+  16  D,  +  D„  .  4-  D .  ,  +  D,  , ,  +  D .  ,  .  (21) 

l  ,m  i,m-\  JL-l,m  &,m  + 1  l+l  ,m  v  ' 


The  system  of  linear  equations  may  be  expressed  in  matrix 
notation: 


M  P 
«  — 


(22) 


The  arbitrary  row  of  the  coefficient  matrix  M  ,  that  row  corresponding  to 
the  H  ,m  element,  may  be  exhibited  in  stencil  form,  acting  on  the  two- 
dimensional  grid  array  of  reference  locations  for  the  pressure  values: 
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G  =  +  A,  p  „  -  B  „  U.  +  B .  ,  u  „  , 

Z,m  Z,m  Z,m  Z  ,m  ^Z  ,m  Z,m-\KZ,m-\ 


Z ,mU Z ,m+C l-l,m  vZ-\,m 


-  4D .  L  „ 
Z,m  Z ,m 


+  D  .  L .  .  +  D  „  -  L ..  .. 

Z,m- 1  Z,m- 1  Z-\ ,m  Z-l  ,m 

+  ^Z,m  +  l  Li(w+ 1  +  tirt^Z+l  ,m 


(23) 


The  matrix  M  is  symmetric  positive  definite. 

'ft 

The  formal  solution  for  p  is  expressed  by  matrix  inversion: 

p*  =  M_1  G  .  (24) 

—  RS  — 

2  2 

However  since  M,  for  a  125x125  grid  array,  has  dimensions  125  xl25  , 
the  inversion  of  M  is  not  a  feasible  calculation.  Recourse  must  be  made 
to  iterative  methods  for  arriving  at  approximate  solutions  for  p*  .  These 
methods  are  referred  to  as  relaxation  techniques. 

The  blending  solution,  as  exhibited  so  far,  includes  no  explicit, 
reference  for  the  determination  of  the  resultant  weight  field,  A  .  However, 
the  fact  that  the  resultant  p*  elements  are  specified  by  the  blending  system 
implies  that  A*  is  also  specified,  and  apparently  this  specification  is 
implicit  in  the  blending  system.  This  conclusion  assumes  that  the  blending 
system  is  conceptually  consistent.  These  matters  require  some  sleuthing, 
schematic  analysis,  and  verifications. 
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2.4.3  Implied  Specification  of  A 


Equation  (24)  states  that,  for  an  arbitrary  SL  ,m  component. 


p*  =  M .  .  G .  +  X  M.  .  ,  G, 

*SL,m  SL ,m  .  p*  m  i,«; i,J  i,j 


(25) 


where  M, 


is  the  diagonal  element,  and  the 


a  n  u.v  v.iw.  uiwii . ,  biiiw  uie  hi «  ,  (i  ,  j  X ,  Tfi )  are 

«  z  W  /  ^  / »»  «  /  ^  /  i  /  j  j 

the  off-diagonal  elements,  of  the  i,m  row  of  the  inverse  matrix,  M 


R3 


By  separating  G.  into  two  parts, 
&9m 


G.  =  A.  p  +  H. 
l,m  SL  ,m  vSL,m  SL  ,m 


(26) 


Eq.  (25)  can  be  expressed  in  the  form 


pX,w  M SL,m;JL,m  A SL,m  pX,m 


+  M .  .  H .  +  X  M.  .  .  G.  , 

SL,m  :SL  ,m  l,m  L  m  t,m;  i,j  i,j 


If  we  have  obtained  the  desired  resultant  which  is  schematically 
expressed  by  Eq.  (17),  then  the  implications  are  that 


ASL,m  pJL,m  Ai,mpSL,m/'ASL,m 


(27) 


(28) 


and 


M  H  +  /  i  M ,  .  ,  G,  , 

^SL.miSL.m^SL.m  irfSL.m 


a  a  ,  * 
ASL,m  p JLlm/AJL,m 


(29) 
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According  to  Eq.  (28)  the  elusive  A  elements  are  given  by  the  diagonal 
elements  of  the  inverse  matrix: 


M 


-1 


(30) 


Inspection  of  the  form  of  Eqs.  (22)  and  (24)  and  of  the  functional  dependencies 
of  the  matrix  elements  reveals  compatability  with  this  conclusion.  However, 
additional  support  is  desirable. 


2.4.4  The  Calculation  of  A* 


For  large  grid  arrays  it  is  generally  not  feasible  to  attempt  the 
inversion  of  the  matrix  M  .  And,  it  is  just  about  as  unrealistic  to  attempt 
a  frontal  attack  on  solving  only  for  the  diagonal  elements  of  the  inverse ,  M~  . 

The  conceptual  relationship  expressed  by  Eq.  (17)  suggests  a  practical 
approach — optimizing  between  accuracy  and  computation  costs .  According 
to  Eq.  (17)  a  finite  change  in  p.  ,  keeping  all  other  forcing  components 

Jit  ,771 

and  all  weights  fixed,  results  in  a  corresponding  finite  change  in  p^  m  as 
expressed  by 


(31) 


The  scheme  is  applied  after  obtaining  a  complete  resultant  distribution 
of  p*  ,  of  adequate  accuracy. 

The  scheme,  implied  by  Eq.  (31) ,  is  applied  for  the  determination 
of  A*  at  one  reference  location  at  a  time.  Each  such  application  involves 
re-solving  for  the  field  of  p*  .  At  first  glance  this  implies  an  inordinate 
amount  of  computing.  However,  the  region  of  influence  of  p.  on  the  p 
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field,  diminishes  with  distance  from  the  l ,m  location.  Hence  p*  need  be 
re-solved,  subject  to  the  specified  change  6p  ,  only  for  a  limited  region 

Aj.fTn 

about  the  i ,m  location. 

In  areas  abundant  with  A  weight,  the  region  may  be  quite  limited. 

In  other  areas,  where  A  weight  is  sparse,  the  region  should  be  more 
extended  to  obtain  comparable  accuracy. 

The  re- solving  of  p*  ,  subject  to  6p.  ,  determines  6p*  .  And, 

z,m  x/t 

according  to  Eq.  (31), 


SL  ,m 


A  6p  / 6p* 

i  ,m  i,m  l,m 


(32) 


A  similar  scheme  can  be  used  to  calculate  the  combined  weight  of 
the  other  parameters  whose  values  are  defined  by  p*  .  For  example. 


H  ,171 


JL,m  + 1 


&  ,m 


(33) 


Consider  u*  :  to  be  the  resultant  according  to 
rSL,m 


■k  *  _  ,  _a  a 

B  U  =  B  U  +BU 
i,m  l,m 


(34) 


Re-solving  for  p*  ,  subject  to  the  single  finite  change  6  fj  ,  determines 

171 

the  corresponding  finite  change  6ji*  .  Accordingly, 

Xj  t17l 


BJl,m  Bi,m  l,m^ ^ l,m 


(35) 
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2.4.5 


Additional  Support  for  the  Blending  Formulation 


Additional  evidence  that  the  desired  blending  is  effected  by  the 
developed  formulation  for  implicit  blending  is  as  follows: 

Numerical  tests  verify  that  this  implicit  blending  produces  the  same  parameter 
values,  with  the  same  weights,  as  are  produced  by  explicit  blending  in 
single-path  arrays. 

Numerical  tests  verify  that  for  any  arrays  the  ambient  contribution,  p^  of 

weight  A^  ^  ,  stemming  from  all  reference  locations  other  than  l,m  ,  is 

independent  of  p„  and  of  its  weight  A. 

* l ,m  i.m 
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2.5 


Additional  Spreading  Parameters 


2.5.1  Expansion  Terms 

Additional  higher-order  difference  parameters  may  be  added  to  the 
blending  scheme  in  order  to  enhance  the  spreading  of  independent  information 
and  in  order  to  effect  additional  controls  on  the  character  of  the  resulting 
analysis . 

Additional  parameters  are  suggested  by  terms  of  the  extrapolation 
expansion: 


P 


p  +  6  ir  •  V  p  +  ~  6  ir  •  V  6  ir  •  V  p  +  .  .  . 


(36) 


where  the  displacement  is 


6  ir  =  ii  fix  +  jj  6y 


•The  terms  may  be  expressed  in  Cartesian  coordinates: 


6  ir  •  Vp  =  6x  +  6y 

dx  oy 


(37) 


6  Ir  •  V  6  ir  •  Vp 


6x2  ^H“+6y2  M  +  26x6y|-£- 


(38) 


This  expansion  suggests  that  the  second  derivatives  and  the  cross  derivative 
are  significant  in  the  spreading  of  information. 
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2.5.2 


Second-Difference  Parameters 


The  second  differences  are  defined  by  the  parameters 


/s 


Z  ,m 


P Z  ,m  +  \  +  Pi,/72-l 


-  2  p 


Z  ,m 


s\ 

L/  =>  D  +  D 

Z  ,m  ^Z+\,m  p Z-\,m 


-  2  p 


(39) 

(40) 


It  is  adequate  for  present  purposes  that  the  weights  of  relevant  available 

ft.  and  ft,  information  be  equal.  Denote  this  weight  by  F.  .  The 
*Z,m  Z,m  M  Z,m 

contribution  to  the  error  functional  is 


“  ■  £{pi.»W. 

Z  ,m  \ 


m+ 1  +  PZ,m- 1  "  2  PZ 


-  ft  v 

,m  nZ,m/ 


+  F 


Z  ,m  (PX+ 1  ,/7z  +  p£-l ,m  2  P4 


-o,  )2( 

,m  Z  ,m  )  | 


(41) 


The  resulting  additions  to  the  system  of  blending  equations  may  be 
expressed  as  contributions  to  Eq.  (23): 


The  contribution  to  S.  .  given  by  Eq.  (21),  is 

Z  ,m 


A  c  =op  4-P  4*F  4-F  4-  F 

j l  ,m  Z  ,m  Z,m- 1  Z,m  + 1  Z-\,m  Z+l,m 


(42). 
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The  contribution  to  the  forcing  element,  G.  ,  is 

a  ,n 

£  ”  2FZ,m  +  F  l  ,m-l  +  +  +  l 


2V  'O'  +  ¥  'O'  +  F  'O' 

i,m  i,m  A-l#w  l-l,m  i+1  ,m  Jt+1  ,m 


(43) 


2.5.3  Cross-Difference  Parameter 

The  cross  difference  is  defined  by  the  parameter 


^ \,m  P4+1.//7j+1  +  P i  ,m  PX+l/772  .  p4,ttc  + 1 


The  weight  is  denoted  by  K.  .  The  contribution  to  the  error  functional  is 

&  ,771 


6E  = 


*  +  * 
p4+l,/rc  +  l  p  4  ,m 


_  *  _  * 

p4+l  ,m  p4,/7!+l 


-U2}  • 


The  contribution  to  S.  is 

4  ,m 


6S„  =  K .  +  K .  .  .  +  K,  +  K.  . 

i  ,m  4  ,m  4-l,m-l  l,m- 1  1-1  ,m 


The  contribution  to  the  stencil  on  p  is 


i 

i 

»— ■ 

«• 

a 

+  K.  . 

4-1,/tz 

-  K„ 

4  ,m 

+  K. 
i,m 

-  K„  1 
i-l,m 

6S. 

4  ,m 

-  K„ 

4  ,m 

-  K.  . 

4-1 ,m-l 

-  K 

4,772  -1 

+  K.  n  . 

4—  1  ,m  —  1 

-  K.  _ 

4-1  ,m-l 

+  K„ 

4,772  - 1 

-  K. 

4  ,/rc  -1 

The  contribution  to  the  forcing  element.  G.  .is 

z  ,m 


6G.  5  K  y  +  K  y.  , 

Z  ,m  Z,m  Z  ,m  Z-l ,m-l  Z-l,m-l 


-  K.  .  y.  -  K.  .  V.  . 

Z,m-1  Z,m-l  Z-\,m  Z-l  ,m 


(47) 


2.5.4  The  Spreading  Parameters 

We  distinguish  between  the  data  parameters  and  the  spreading 
parameters.  The  data  parameters  are  the  pressure,  p  ,  and  the  first 
differences,  /I  and  v  .  The  spreading  parameters  are  of  higher  order  and 
include  the  Laplacian,  L  ,  the  second  differences,  ^  and  D* ,  and  the 
cross  difference,  7. 

The  information  in  the  spreading  parameters  may  be  considered  to 

/ 

include  information  derived  from  differences  of  these  differences,  that  is, 
from  terms  of  their  analytic  expansions — from  higher-order  terms  to  which 
there  is  no  explicit  reference  because  of  truncation.  Thus  their  weights 
do  not  represent  independent  information:  the  variances  are  spatially 
correlated.  However  this  does  not  violate  the  blending  conceptions.  The 
information  in  the  spreading  parameters  may  be  considered  to  be  the 
resultant  of  a  blending  operation  which  assimilates  all  higher-order 
information,  preparatory  to  the  next  step  in  downward  blending  which  is 
actually  performed.  This  conception  defines  blending  to  be  an  integrating 
procedure  which  may  be  performed  stepwise  in  order. 

The  spreading  values  and  weights  have  spatial  continuity.  In  a 
sense  they  are  analyzed  fields.  The  local  addition  of  information  to  a 
spreading  parameter  should  be  effected  by  a  reanalysis,  or  restructuring, 
procedure  based  on  the  analytical  continuity. 
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2.6 


The  general  method  of  Successive  Over-Relaxation  (SOR)  is 
appropriate  for  the  solution  of  the  linear  system  expressed  by  Eq.  (22) . 
The  matrix  M  is  normalized;  this  procedure  makes  each  diagonal  element 
unity.  The  matrix  is  factored: 


M  =  I  +  U  +  V  (48) 

ft*  &  &  & 


where  I  is  the  identity  matrix,  U  has  non-zero  elements  only  above  the 
main  diagonal,  and  V  has  non-zero  elements  only  below  the  main  diagonal. 

ft* 

Equation  (22)  expands  into 

p*  =  G  -  Up*  -  Vp*  .  (49) 

—  ~  —  pa  ~ 

The  direct  SOR  method  is  defined  by 

u  p*(r) 

ft*  — 

where  the  superscript  in  parentheses  denotes  the  pass  (i.e.,  iteration) 
number  of  the  successive  estimate.  The  method  converges  for  0  <  uj  <  2. 
Over-relaxation  implies  0)  >  1 .  We  are  still  free  to  design  the  ordering 
of  the  elements  of  the  vector  p*. 

Of  the  orderings  tested,  the  following  ordering  performed  best — 
that  is,  it  appeared  to  give  the  most  rapid  convergence. 

The  elements  of  p*  are  divided  into  five  subsets.  A  subset 
includes  every  fifth  element  of  each  row  of  the  grid-point  array,  with  a 
displacement  of  two  elements  from  one  row  to  the  next.  The  division  is 
illustrated  by  Fig.  7. 


*  (r+ 1)  _  *(r) 


+  (A)  (  G 


-  f*(r)) 
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Fig.  7  Labelling,  as  to  subset,  of  a  section  of  the  grid-point  array. 


The  ordering  of  the  vector  p*  is  specified  as  a  sequencing  of 
the  subsets: 


/  1  \ 

2 

£ 

4 

\~/ 


(51) 


The  significance  of  this  separation  is  that  the  ordering  within  each 
subset  is  immaterial.  The  calculation  of  an  element  is  independent  of  the 
other  elements  of  the  same  subset  as  that  element.  In  effect,  in  computing 
each  subset,  calculated  according  toEq.  (50),  the  scheme  is  effectively 
both  successive  and  simultaneous.  However,  in  a  full  pass  encompassing 
all  five  subsets,  the  scheme  is  definitely  successive. 

The  first  full  pass  is  calculated  with  U)  =  1.  In  the  second  pass 
U)  is  set  high,  »  1 . 95 ,  and  is  decreased  with  each  successive  pass .  The 
last  pass  is  made  with  set  at  1. 
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2.6.2  Calculation  of  A* 

•jg 

For  each  element,  the  calculation  of  A.  involves  a  re-solution 
•  Ji  ,m 

of  p*  subject  to  an  adjustment  in  the  forcing  element,  6p.  .  This 

J i  ,m 

re-solution  is  performed  only  for  a  limited  region  about  the  l ,m  location. 

In  areas  abundant  with  A  weight,  the  region  may  be  quite  limited.  In 
other  areas,  where  A  weight  is  sparse,  the  region  should  be  more  extended 
to  obtain  comparable  accuracy. 

The  method  adopted  for  this  re- solution  of  p*  is  the  SOR  method 
but  with  U)  =  1, — that  is,  we  use  Eq.  (49): 


p*(r+l)  _  q  _  u  p*(r)  _  v  p*(r+l) 

~  ~  RS  ~  RS  ~ 


(52) 


The  ordering  for  the  limited  region  is  by  sequencing  as  to  subset,  with 
subset  identification  according  to  Fig.  7. 

Before  beginning  with  the  iterations,  using  Eq.  (52),  the  first-guess 
for  p*  is  adjusted  at  all  elements  of  the  limited  region,  by  the  increment 


6p* 


1 


where  A  is  a  specified  approximation  of  background  weight. 

B 

The  elements  of  p*  which  surround  the  limited  region  are  not 

changed  in  response  to  6p.  .  This  is  a  restraint  on  6p.  and  causes 

i i,m  SL  ,m 

an  overestimate  of  A*  .  The  first-guess  adjustment,  according  to 

& ,  m 

Eq.  (53),  introduces  a  compensating  bias  which  is,  however,  dissipated 
if  the  limited  region  is  taken  very  close  to  convergence. 
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2.7 


The  Time  Axis 


2.7.1  Explicit  Use  of  First-Guess  Fields 

How  should  the  value  inherent  in  a  first-guess  field,  predicted 
from  the  preceding  analysis,  be  exploited?  Wherein  lies  the  information? 

In  the  gradients?  In  the  Laplacian?  For  example,  a  500-mb  barotropic 
model  may  be  quite  skillful  in  predicting  the  Laplacian  distribution,  but 
may  have  considerably  more  variance  in  predicting  the  absolute  height 
distribution.  How  should  each  of  the  parameter  information  fields  be 
weighted  ? 

It  is  easier  to  conceive  of  first-guess  fields,  for  component 
parameters,  to  be  weighted  as  to  their  resolution  rather  than  as  to 
independent  worth.  The  resultant  resolutions  realized  for  the  preceding 
analysis  decrease  with  increasing  variance  growth  due  to  error  drift  in 
prediction  skills. 

Resultant  weights  are  appropriate  for  spreading  parameters.  However, 
first-guess  fields  for  data  parameters  must  be  weighted  according  to 
independent  worth  for  assimilation  with  the  new  observations. 

If  all  parameters  are  derived  from  a  single  resultant  prediction  of 
the  object  parameter,  then  all  the  parameters,  p  ,  \i  ,  v  ,  L  ,  etc. ,  are 
mutually  consistent.  However,  such  consistency  is  not  essential. 

The  problem  of  determining  the  independent  worth  of  first-guess 
fields  for  data  parameters  reduces  to  the  following: 

Given  A*  ,  B*  ,  C*  ,  D*  ,  determine  A  ,  B  ,  C  . 
o  o  o  o  '  o  o  o 

This  inversion  conceivably  gives  rise  to  negative  elements  in  the  independent 
weights.  This  negative  extension  of  information  weight  is  formally  tractable 
and  presents  no  difficulties.  We  have  not  been  able  to  solve  this  inversion 
problem  except  for  trivial  examples. 
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The  procedure  assumed  for  initial  weighting  of  first-guess  fields 
is  assignment  by  specially  designed  formulas. 

While  expedient,  this  direct  approach  is  a  far  cry  from  the 
comprehensive  use  which  could  be  made  of  data  along  the  time  axis,  if  it 
were  not  for  computer  capacity  limitations . 


2.7.2  Comprehensive  Approach 
2.7.2. 1  Explicit  Use  of  Predicted  Change 

The  comprehensive  approach  lies  in  the  design  of  an  error  functional 
which  couples  several  analysis  times,  e.g., 

t  =  ...  (n-2)r  ,  (n-l)r  ,  nT  (54) 

where  T  is  the  period  between  analysis  times.  The  blending  is  thereby 
extended  to  three  dimensions. 

At  the  new  analysis  time,  nr  ,  the  error  functional  may  couple  the 
times  nr  ,  (n-l)r  and  (n-2)r  .  The  information  which  couples  these 
analysis  times  is  in  the  change  fields,  6p  ,  6/lx  ,  ,  etc. ,  obtained  from 

a  prediction  scheme  calculating  on  earlier  analyses  from  time  (n-2)r  to 
(n-l)r  and  from  time  (n-l)r  to  nr  .  These  predicted  changes  are 
appropriately  weighted  as  data  parameters  and  as  spreading  parameters. 

The  three-dimensional  blending  produces  resultant  analyses  for 
times  (n-2)r  ,  (n-l)r  and  nr  .  In  an  operational  configuration  the  resultant 
for  the  latest  time,  nr  ,  is  the  new  current  analysis.  The  resultant  for  the 
preceding  time,  (n-l)r  ,  is  an  update  analysis  which  includes  off-time  data 
both  ahead  and  behind — achieving  forward  and  backward  continuity  in  time. 
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2. 7. 2. 2  Implicit  Dynamics 


The  comprehensiveness  of  the  analysis  system  would  reach  ultimate 
levels  by  making  the  predicted  changes,  6p  ,  6/j, ,  ,  etc. ,  implicit 

functions  of  the  analyses,  running  the  governing  dynamics  both  forward 
and  backward  in  time.  The  resultant  analyses  would  be  free  of  the  problem 
of  update  shock,  having  achieved  optimum  adjustments  of  continuity  in 
time. 
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3. 


Additional  Details  of  FIB/SLP 


In  this  section  further  details  of  the  FIB/SLP  operations  are 
discussed.  A  basic  flow  chart  is  shown  in  Fig.  8. 

3 . 1  First-Guess  Field  Preparation 

The  first  guess  to  the  hemispheric  sea-level  pressure  analysis  is  a 
combination  of  the  applicable  prognosis  field  and  the  previous  hemispheric 
and  tropical  analyses  extrapolated  by  a  kinematical  technique. 


3.1.1  Combination  of  Previous  Analyses 

Grid  point  pressure  values  are  first  interpolated  from  the  previous 
tropical  analysis  at  each  hemispheric  grid  point  within  the  tropical  grid 
(i .  e . ,  south  of  60°  N  latitude) .  A  weighted  combination  of  the  tropical 
interpolated  and  hemispheric  values  is  made  by  the  formula 


Po  =  wPt  +  (1— w)  ph 


(55) 


where  w  =  0.86  -  sin  <p  ,  (0  §  w  =  1)  , 

<p  is  latitude, 

p,  is  the  hemispheric  grid  point  pressure  value, 
h 

Pt  is  the  interpolated  tropical  value,  and 
Pq  is  the  resulting  weighted  value. 

The  formula  gives  a  weight  to  the  tropical  value  of  zero  at  60° N,  0.36  at  30° N, 
0.86  at  the  equator  and  1.0  from  8°S  to  the  corners  of  the  grid  (19° S). 

Use  of  the  tropical  analysis  provides  continuity  with  southern 
hemisphere  data  and,  implicitly,  with  climatological  information. 


-43- 


Fig.  8  Schematic  Flow  of  Operations  for  FIB/SLP 
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*  £ 


3.1.2 


Kinematical  Extrapolation 


Rather  than  utilize  the  six-hour-old  combined  analysis  directly,  an 
extrapolation  is  applied  to  minimize  the  effects  of  system  movement.  The 
residual  500-mb  height  analysis  (SR500) — the  actual  500-mb  height  analysis 
with  the  disturbance  scale  removed — is  used  as  a  steering  field.  From  the 
SR500  height  values  geostrophic  winds  are  used  to  derive  six-hour  movements. 
The  program  then  looks  upstream  from  each  grid  point  in  the  six-hour-old 
combined  surface  analysis .  The  distance  is  specified  by  a  percentage  of 
the  six-hour  movement.  The  pressure  interpolated  there  is  brought  back  to 
the  grid  point.  The  optimum  speed  has  been  found  to  be  about  58%  of  that 
specified  by  the  SR500  geostrophic  flow  with  a  daily  variation  of  about  10%. 

A  continuing  evaluation  of  the  optimum  percentage  is  made  by  the  program. 


3.1.3  Prognostic  Field 

The  FNWC  primitive-equation  surface  pressure  forecast  is  used  as 
the  prognostic  portion  of  the  first  guess.  Depending  on  the  analysis  time, 
this  is  either  a  6-hour  or  12-hour  forecast,  while  the  kinematical  extrap¬ 
olation  is  always  for  a  6-hour  period.  These  two  fields  are  combined  to 
produce  the  first  guess  (pq)  to  FIB/SLP.  The  prognostic  portion  is  weighted 
as  a  function  of  the  sine  of  latitude  so  that  in  southern  latitudes,  the 
extrapolated  field  carries  relatively  more  weight.  The  total  weight  of  this 

field  (A  )  increases  from  0.001  at  the  north  pole  to  0.01  south  of  the 
o 

equator.  To  account  for  the  information  in  the  tropical  analysis  outside 
the  borders  of  the  hemispheric  grid,  the  weight  values  are  increased  on 
and  one  row  inside  the  border. 
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3.1.4 


Higher  Derivative  Fields 


From  the  first-guess  pressure  field  computed  above,  the  first-guess 

finite  difference  (u  ,  V  ) ,  Laplacian  (L  ),  second-difference  ('tt  ,  'O'  )  and 

o  o  o  o  o 

cross-difference  (7  )  fields  are  calculated  from  the  formulas  in  Fig.  2  and 

Eqs.  (39),  (40)  and  (44).  The  Laplacian  weight,  D  ,  is  computed  as 

o 


o  a  ,m 


f  +  f,  L  * 
o  1  o  SL,m 


(56) 


2  2 

where  L  „  is  the  mean  of  L  at  grid  point  ,m)  and  the  four 
O  *  ,772  O 

surrounding  grid  points  and  f  and  f,  are  adjustable  constants  (Appendix  1) . 

o  l 

The  first-difference  weights  (B  ,  C  ) ,  second-difference  weights  (F  )  and 

o  o  o 

cross-difference  weights  (K  )  are  computed  as  constant  factors  of  D  .  All 

o  o 

weights  are  set  to  zero  where  their  influence  extends  outside  the  grid,  e.g. , 
Dq  at  all  boundary  points,  Bq  along  the  top  boundary,  and  Cq  along  the 
right  side  of  the  grid.  (Bq  =  Cq  at  all  interior  grid  points.)  The  Bq  and  Cq 
values  are  increased  where  their  respective  /i  and  v  values  refer  to  a 
boundary  point;  this  extends  the  shape  of  the  tropical  field  into  the  boundary 
regions  of  the  hemispheric  analysis. 


3 . 2  Assembly  of  New  Information 

New  pressure  and  wind  data  are  assembled  from  ship  (SH) ,  land 
six-hourly  synoptic  (SM)  and  land  hourly  (SA)  reports.  The  data  lists  are 
merged  and  sorted  so  that  an  hourly  report  will  be  used  only  if  there  is  no 
synoptic  report  available  at  the  same  location.  When  the  report  at  analysis 
time  is  missing,  the  one  closest  to  the  analysis  time  (up  to  +  3  hours)  is 
admitted  to  the  data  list  along  with  an  indication  of  report  age.  This  may 
bring  in  three-hourly  (SI)  reports,  but  at  considerably  reduced  weight,  as 
discussed  later. 
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3.2.1 


Pressure  Reports 


The  pressure  data  are  reported  in  shortened  form  with  the  100's  of 
millibars  omitted.  The  first  step  is  to  determine  the  multiple  of  lOOmb 
that,  when  added  to  the  report,  will  minimize  the  difference  between  the 
first  guess  and  the  modified  report.  As  a  result  of  this  operation,  no 
report  will  differ  from  the  first  guess  by  more  than  50 mb. 

The  pressure  report  can  now  be  extrapolated  to  the  nearest  grid 
point,  determined  by  rounding  the  exact  i,j  coordinates  to  the  nearest 
integers .  Extrapolation  is  made  along  the  gradients  of  the  best  available 
field,  i.e. ,  the  pQ  field  if  the  analysis  is  the  first  one  being  performed  in 
the  synoptic  period,  the  result  of  first  pass  blending  if  in  the  second  pass, 
or  the  final  analysis  for  the  same  time  period  if  performing  an  update 
analysis.  The  formula  for  extrapolating  pressure  reports  assumes  that  the 
difference  between  the  first  guess  and  the  report  remains  constant  over 
the  distance  from  the  report  to  the  grid  point: 

p H  ,m  (Po)j i  ,m  +  (Pn  Po)i,j 

where  p„  is  the  extrapolated  pressure  at  grid  point  (X  ,m) ,  (p  )  .  _  is 
sl  ,m  O  X,  ,171 

the  guess  value  at  the  grid  point,  (p  ,  is  the  reported  value  at  station 

n  i, j 

location  (i,j),  and  (pQ).  j  is  the  guess  value  interpolated  at  the  station 
location  using  a  Bessel  16-point  operator. 

The  maximum  allowable  difference  between  the  report  and  best 
available  guess  increases  with  the  magnitude  of  the  gradients  in  the 
vicinity  of  the  report  and  with  latitude.  The  report  is  considered  a  gross 
error  and  rejected  if 

lpn-poU.M  >  f2  +  f3  +  f4  »  '  cos  -  <S8> 


-47- 


where  |  g  |  ^  is  the  mean  absolute  gradient  in  the  surrounding  grid  points, 
<p  is  latitude  60°)  and  f ^  and  f  are  adjustable  constants.  In  an 

area  of  strong  gradients,  north  of  60° N  latitude,  the  allowable  tolerance  may 
reach  30mb. 

The  reliability  of  a  pressure  report  decreases  with  the  magnitude  of 
gradients  in  the  vicinity  of  the  report,  with  elevation,  and  with  age.  In 
terms  of  variance 


2 

CT 

n 


+  fc  g. 

6  a  ,m 


2 

+  f?  E  + 


.f8  +  f9  gj2,w- 


(59) 


where  a  is  the  total  variance  of  the  pressure  report,  g.  is  the  mean 
n  J l  ,m 

square  gradient  surrounding  the  report,  E  is  elevation,  H  is  age  and  f 

o 

through  fg  are  adjustable  constants.  The  class  of  report  is  represented 

in  f  .  The  age  term  does  not  affect  the  intrinsic  variance  of  a  report,  but 
o 

does  reduce  its  weight  for  assembly  purposes.  The  weight,  or  reliability, 

-2 

of  the  report  is  A  =  ct  where  n  refers  to  a  single  (unassembled) 
n  n 

reliability.  If  several  conflicting  reports  are  received  from  the  same 
station,  their  individual  weights  must  be  reduced  since  they  are  not 
independent  information — their  total  weight  should  be  . 

The  assembled  value  of  pressure  at  a  grid  point  is  a  weighted  mean 
of  all  data  referred  to  the  grid  point,  with  each  value  contributing  to  the 
mean  in  accordance  with  its  reliability.  Such  a  weighted  mean  can  be 
derived  one  report  at  a  time,  independent  of  assembly  order,  by 


=  (p„  ,Am  ,  +  p  A  )/( 

\  N-l  N-l  n  n/  V 


Am  1  +  A  ) 

N-l  n/ 


where  p  ,  A  refer  to  a  single  report  and  P._  ,  A.,  refer  to  the  assembly 
n  n  N  N 

of  N  reports. 


-48- 


The  assembly  is  initialized  with  the  first-guess  pressure  field  •  . 

(p  ,  A  ) .  This  must  be  distinguished  from  the  best  available  field — which, 
o  o 

at  least  partly,  represents  current  information  and  cannot  be  considered 
independent.  An  example,  considered  later,. is  an  analysis  made  earlier 
for  the  same  time  period.  The  first-guess  pressure  field  is  derived  totally 
from  past  information,  extrapolated  forward  in  time,  and  is  therefore 
independent  information ,  validating  the  addition  of  variance  rule  (Section  2.1.2). 


3.2.2  Wind  Reports 

In  order  to  use  a  report  of  wind  direction  and  speed,  conversion 
must  first  be  made  to  components  of  pressure  finite  differences  along  the 
grid  coordinates  by  an  appropriate  balance  equation.  The  balance  equation 
is  written  in  vector  form  as 

W  .  V  W  +  f  Ik  x  W  +  X  W  =  f  Ik  x  W  ,  (62) 

where  W  .  is  the  geostrophic  wind,  W  is  the  actual  wind,  Ik  is  the  unit 
9 

vertical  vector,  x  is  a  frictional  parameter,  and  f  is  the  coriolis  parameter. 
Local  accelerations  of  the  wind  are  omitted.  The  frictional  effect  on  the 
surface  wind  is  contained  in  the  x  term.  In  component  form  the  geostrophic 
wind  may  be  written  as 


(63) 
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where  the  derivative  subscripts  refer  to  geostrophic  wind  components  in 
the  best  available  field  (e.g.,  an  analysis  for  the  same  time)  and  u,v  are 
the  reported  wind  in  component  form.  The  value  of  R  ,  slightly  greater 
than  1.0,  has  been  added  to  increase  the  difference  between  the  geostrophic 
and  actual  velocities  without  increasing  the  inflow  angle.  This  empirically 
accounts  for  the  fact  that,  due  to  vertical  momentum  transfer,  the  frictional 
force  is  not  directed  exactly  opposite  the  surface  wind  but  rather  opposite  a 
direction  somewhere  between  the  surface  and  free  flow. 

The  x  term  allows  for  variations  in  roughness.  It  can  be  shown 
from  Eq.  (63),  assuming  negligible  curvature,  that  x  =  f  •  tan  a.  where  a 
is  the  inflow  angle  of  the  surface  wind.  The  value  of  X  is  set  for  an 
approximate  inflow  angle  of  12°  over  the  sea  and  22°  over  land. 

The  derivative  terms  are  computed  as  four  separate  fields  using  a 
16-^point  operator.*  The  values  of  the  derivatives  are  then  determined  for 
each  report  by  interpolation  to  exact  station  location. 

The  final  pressure  finite  difference  values,  and  ,  are 
computed  by  applying  the  geostrophic  wind  equation  to  the  u  and  v 
components  found  above,  and  solving  for  pressure  differences  over  1  grid 
length  in  the  m  and  l  grid  directions ,  respectively. 

The  variance  of  the  wind  report  is 


n 


=  f10  sin2? 


1  +  sin  (/? 

1.866 


,  2 
f6  qi,m 


f7E 


f8  +  f9  qi,m 


H 


(64) 


*Holl,  M.  M. ,  1966;  "The  assimilation  of  wind  observations  in  the  analysis 
of  the  atmospheric  mass  structure"  ,  Quarterly  Report  1,  Contract  No. 
N00228-66-C-1923,  (Fleet  Numerical  Weather  Central) ,  Meteorology 
International  Incorporated,  Monterey,  California,  16  pp. 
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where  f ^  is  a  class  variance.  The  first  term  on  the  right  is  the  "calm 

wind"  variance.  The  sin^<p  portion  is  needed  to  convert  from  a  wind 

variance  to  a  variance  in  terms  of  pressure  difference  and  the  remainder 

of  the  term  accounts  for  variation  in  grid  spacing  with  latitude.  Note  that 

a  light  wind  in  the  tropics  can  have  a  very  small  variance  when  expressed 

in  terms  of  an  equivalent  pressure  difference.  Since  pi  and  V  are 

n  n 

referred  to  the  same  grid  point,  the  same  weight  can  be  assigned  to  each, 

_2 

i.e.,  B  =  C  =a  .  Separate  fields  of  B  and  C  are  maintained  because 
n  n  n 

of  different  boundary  values  on  the  sides  where  pi  or  V  would  extend 
outside  the  grid. 

The  u  and  v  pressure  differences  and  their  weight  values 
n  n 

(B  =  C  )  are  assembled  in  the  same  manner  as  the  pressure  reports, 
n  n 

beginning  with  first-guess  fields  piQ  ,  v  ,  Bq  and  Cq  .  No  extrapolation 
is  performed.  The  reference  grid  point  which  minimizes  the  pi  and  v 
position  error  is  found  by  subtracting  0.25  from  the  exact  i  and  j  coordinate 
and  rounding  the  result.  This  offset  is  necessary  when  referring  both 
differences  to  the  same  grid  point  considering  the  different  locations  of  the 
fi  and  v  reference  points  (See  Fig.  2). 

No  initial  gross-error  rejection  is  performed  other  than  a  limit  on 
reported  wind  speed  of  60m/sec. 

Upon  completion  of  the  pi  and  v  assembly,  the  assembled  weights 
are  reduced  at  grid  points  containing  data  to  account  for  the  uncertainty  in 
the  balance  equation  conversion.  (See  Section  2. 3. 2. 3.)  Expressed  as  a 
variance: 


(65) 


where  o 
a 

of  B  or  C) 


is  the  assembled  variance  at  the  grid  point  (i.e. ,  the  inverse 
o 

,  a  is  the  variance  of  the  balance  approximation,  and  cr 
b  bo 
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is  the  final  variance  at  the  grid  point  from  which  B  and  C  are  derived. 
Thus ,  the  total  weight  at  a  grid  point  is  bounded  by  the  balance 
approximation  variance,  no  matter  how  many  wind  reports  are  assembled  at 
the  grid  point.  This  also  reduces  weights  at  grid  points  with  single  wind 
reports.  The  assembled  reliabilities  included  variance  in  the  wind  report 
only.  In  terms  of  pressure  difference,  this  variance  may  have  been  quite 
small,  resulting  in  large  assembled  weight  before  reduction  for  the  balance 
approximation . 


3 . 3  Blending  for  p* 

When  all  data  have  been  assembled  into  their  respective  fields,  the 
program  can  blend  the  information  into  a  best-fit  p*  field  which  minimizes 
the  error  functional  (Eqs.  18,  41,  and  45).  The  requirement  of  minimizing 
the  errors  between  the  pressure  field  and  its  derivative  fields  at  all  grid 
points  leads  to  the  basic  implicit  blending  equation  (Eq.  20) .  The  implicit 
nature  of  the  problem  requires  an  iterative  solution.  Considerable 
computational  simplification  is  possible  by  combining  terms  which  are 
constant  throughout  the  iteration  or  which  are  common  multiples  of  the 
same  p*  .  This  leads  to  the  forcing  term  (Eq.  23)  which  involves  the 
starting  fields  only,  and  a  sum  term  (Eq.  21)  which  concerns  only  p*  at 
the  grid  point  being  solved.  Other  terms  involve  p*  at  surrounding  grid 
points  and  can  be  grouped  into  a  stencil  operating  on  the  p*  field.  The 
inclusion  of  second-difference  and  cross-difference  fields  requires 
additional  terms,  as  shown  in  Section  2.5.2  and  2.5.3.  The  total  stencil 
is  shown  below: 
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By  transforming  into  the  following  fields: 

X  =  C .  +  4D  +  4D  +  2F  +2F  +  K  +  K  ,  , 

4  ,m  i,m  JL,m  JL+l,m  4  ,m  4+1  ,m  4,/rc-l  JL,m 

Y.  =  B.  +  4D  +  4D  +  2F  +  2F .  +  K.  +  K  , 

4fw  4,/?j  4,/?j  £,m-l  4,w+l  x>-1  ,m  &  ,m 


Z.  =  D  +  D  +  K  , 

4,/fl  4,/fl  +  l  4+l,/?l  4,/fl 
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®A,m  D i,m  +  F4,/tz  ' 

R4,/rc  ~  D4,/rc  +  D4+l,m+l  +  K4,«  ' 

the  stencil  reduces  to 
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The  sum  term  is: 


S.  =  A  +  B.  +  B.  ,+C,  +  C,  , 

4,772  A  ,m  A,m  A,m- 1  A  ,m  4-1,772 


+  16  D  +  D  .  +  D  ,  +  D  .  f .  +  D  . , . 

A  ,m  A,m- 1  4-1,772  A,m  + 1  4+1  ,m 


+  8  F.  +  F  .  .  +  F  .  +  F .  +  F .  , 

4,772  4  ,772  - 1  4-1  ,m  A,m  + 1  4+1  ,m 


+  K.  +  K.  .  .  +  K.  ,  +  K„  _ 

1,171  4-  1,772-1  A,m- 1  4-1,772 


(66) 


and  the  forcing  term  is: 


G4,772  A4,772  P4,772  ~  BA,m^  A  ,m  +  B4  ,772  - 1  ^A,m- 1  "  C  A  ,m  UA,m  +  C4-l  ,m  U  A-l,m 

-4D  L  +D  .  L.  .  +  D .  .  L.t  +  D  L .  +  D  L ..  , 

A,m  A, m  A,m-l  A,m-l  4-1 ,m  4-1,772  A,m  + 1  A,m+ 1  4+1 ,m  A+l ,m 

-  2F  'Ll'  +  F  ft  +F  /l£ 

A ,m  r A ,m  A,m- 1  *A,m- 1  A,m+ 1  H‘A,m  + 1 

—  ?F  O  j.  p  O  j.  p  O 

A ,m  A ,m  A-l ,m  A-l,m  A+l,m  A+l ,m 

+  k.  y.  +  Kfl  .  .  y.  .  _  -  k„  .  y,  _  -  k.  .  y,  .  ,  (67) 

4,772  4  ,772  4-  1,772  -  1  4-  1,772  -  1  4  ,772  -  1  4  ,772  -  1  4-  1  ,772  4-  1  ,772 

where  '/£'  and  are  the  second  differences  in  the  m  and  4  directions, 
respectively  (Eqs.  (39)  and  (40))  and  7  is  the  cross-difference  (Eq.  (44)). 

Thus,  the  blending  operation  requires  eight  125x125  fields:  X  ,  Y  ,  Z  , 

Q  ,  R  ,  S  ,  G  and  p*  . 
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Because  of  central  memory  storage  limitations,  these  fields  are 
maintained  in  extended  core  storage  (ECS)  and  brought  into  central  memory 
nine  rows  (125  points/row)  at  a  time,  throughout  the  iteration.  (Only  p* 
need  be  stored  in  ECS  after  completing  a  set  of  rows.)  Since  grid-point 
values  two  rows  away  from  the  computation  row  are  required ,  the  nine 
rows  include  four  rows  of  overlap.  This  yields  five  rows  for  actual 
computation,  or  25  sets  to  complete  the  grid  array.  The  p*  computation  is 
made  only  every  fifth  point  across  a  row  in  order  to  realize  the  maximum 
convergence  rate  (See  Section  2.6.1).  The  starting  point  in  each  row 
follows  the  cycle  (1,  3,  5,  2,  4).  Five  passes  through  the  field  are 
needed  before  p*  at  every  point  has  been  modified  once. 

After  the  first  five  passes  (representing  one  full  pass) ,  overrelaxation 
is  applied  to  further  improve  the  convergance  rate.  This  is  done  by  over¬ 
correcting  p*  ,  that  is ,  by  changing  p*  at  a  grid  point  by  the  computed 
change  times  an  overrelaxation  factor.  This  factor  must  be  less  than  2.0 
to  maintain  computational  stability.  The  factor  is  decreased  from  1.95  at 
pass  6  to  1.0  in  the  final  pass.  This  variation  gives  an  even  convergence 
rate  over  the  spectrum  of  wavelengths  encountered  in  the  pressure  field. 

The  first  guess  to  the  p*  iteration  maybe  the  assembled  p^  field 
or  another  better  field,  if  available.  For  example,  when  performing  an 
update  analysis,  the  previously  blended  analysis  for  the  same  time  period 
may  be  used,  which  further  enhances  convergence.  Note,  however,  that 
to  avoid  a  compounding  of  information,  such  a  field  may  be  used  only  as  a 
first  guess  to  the  iteration  and  may  not  be  used  in  the  forcing  term. 

Utilizing  the  methods  discussed  above  to  enhance  the  convergence 
rate,  a  total  of  80  passes,  operating  on  every  fifth  point,  is  found  to  be 
adequate.  This  is  approximately  equivalent,  in  computing  time,  to  16 
passes  when  operating  on  every  point.  A  well-converged  p  field  is 
required  for  the  next  step — computation  of  the  grid-point  resultant  reliability 
field ,  A*  . 
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3.4  Computation  of  A* 


The  pressure  reliabilities  (A)  assembled  at  each  grid  point  represent 
only  the  reliability  of  the  reports  if  they  had  no  effect  on  each  other,  that 
is,  if  there  were  no  blending  or  interaction  between  the  grid  points.  The 
actual  reliability  of  blended  grid-point  pressure  values  should  always  be 
greater  than  the  assembled  reliability  to  account  for  information  propagated 
from  surrounding  grid  points  through  gradient  and  higher-order  information. 
For  example,  the  reliability  of  a  blended  pressure  value  at  a  grid  point  with 
no  assembled  reports  should  be  greater  than  the  first-guess  reliability  if 
reports  surround  the  grid  point  and  if  there  is  gradient  knowledge. 

A  perturbation  technique  is  used  to  compute  A*  by  measuring  the 
leverage  of  a  pressure  change  on  the  information  at  a  grid  point  and  its 
surrounding  area.  To  obtain  an  exact  solution  of  A*  it  would  be  necessary 
at  each  grid  point  to  re-solve  for  p*  over  the  entire  field,  given  a  fixed 
change  of  p  in  the  forcing  term.  The  value  of  A*  would  then  be 


A 


* 

&  ,m 


(68) 


as  shown  in  Eq.  (31).  Stated  in  words,  if  the  assembled  pressure  change 
(6p)  has  small  leverage  and  changes  p*  only  slightly  then  6p*  will  be 
small  and  A*  will  be  large. 

To  keep  the  program  computation  time  within  reasonable  limits, 
the  re-solution  is  carried  out  only  in  the  vicinity  of  the  grid  point  {&,m) 
as  shown  in  Fig.  9. 
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Fig.  9  Computational  sub-area  for  A  indicated  by  (°)  grid  points 

J </  ,171 

Changes  are  made  to  p*  only  at  the  (o)  points.  At,  and  outside, 

k 

the  (x)  points,  no  computation  is  performed  when  solving  for  A  at  grid 
point  (j i ,m) .  Holding  p*  fixed  at  the  (x)  points  is  the  same  as 
assuming  infinite  reliability  at  these  points.  The  effect  is  to  limit, 
somewhat,  the  change  in  p*  and  therefore  overestimate  A*.  This 

Xj  ,171 

overestimation  is  noticeable  only  in  data- sparse  areas  or  when  no 
information  other  than  first-guess  is  found  within  the  (x)  points.  In 
these  situations,  the  computational  area  could  be  enlarged  to  include  at 
least  some  minimum  amount  of  assembled  information.  Computation  of 

k 

A  may  extend  to  the  grid  boundaries.  Points  in  the  sub-area  which  are 
outside  the  grid  do  not  affect  the  result  since  all  weights  extending 
across  the  boundaries  have  a  value  of  zero. 

The  solution  for  A*  is  found  by  applying  the  p*  formula  to  the 
fields  used  in  the  p*  computation  except  that  p*  is  now  a  blended  field . 
The  importance  of  a  well-converged  p  field  is  seen  by  considering  the 
A*  formula — the  change  in  p*  should  result  only  from  the  perturbation 
in  p  and  not  reflect  further  convergence  in  the  original  p* . 


The  nine-row  portions  of  the  fields  are  sufficient  for  computing  A* 
along  only  one  row,  requiring  125  field  exchanges.  At  each  grid  point 
{&  ,m) ,  p  is  changed  in  the  forcing  term  only  at  (&  ,m) .  Then  p*  at  all 
computation  points  is  increased  toward  the  expected  result  to  speed 
convergence .  The  iteration  for  the  modified  p*  value  is  made  through  the 
area  in  a  manner  analagous  to  that  used  in  the  p*  blending. 


©  5  3  1  © 
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Fig.  10  Subset  labelling  in  the  A  computational  sub-area .  Grid  point 
(4. ,m )  is  underlined.  ,m 

In  Fig.  10,  the  #1  subset  points  are  operated  on  first,  then  the  #2 
points,  etc. ,  through  the  #5  points.  (The  circled  points  are  shown  to 
complete  the  pattern  only;  no  computations  are  made  there.)  This  same 
pattern  is  repeated  4  times.  As  shown  in  Section  2.6.1,  ordering  within 
a  subset  is  immaterial.  In  the  last  pass,  the  #5  point  in  the  top  row  is 
skipped  and  a  final  computation  is  made  at  the  center  grid  point  {l,m) , 
underlined.  (Note  that  the  point  skipped  is  not  contained  in  the  stencil 

for  the  center  grid  point.)  A  may  now  be  computed,  knowing  the 

A*  .  171 

#  "k 

6p.  applied  to  the  forcing  term,  the  A  value  and  6p  computed 
SL ,  m  *  ,m  z  ,Th 

aS  ^modified  ^  ^ & ,m  ' 
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Before  proceeding  to  the  next  grid  point,  the  forcing  term  at  {i,m  ) 
is  restored  to  its  original  value  and  A*  is  packed  with  another  field.  The 
p  values  do  not  have  to  be  restored  since  they  are  extracted  for  each  point 
from  the  full  p*  field  which  is  never  modified.  Upon  completion  of  a  row, 
the  packed  A  values  are  written  in  ECS. 

No  increase  in  the  A*  computation  area  is  attempted  in  data 
sparse  regions.  It  has  been  found  that  far-removed  data  have  only  a  small 
effect  on  A  which  is  nearly  masked  by  the  background  reliability.  This 
background  reliability  results  from  the  higher-order  derivatives  which  have 
small  weight  at  individual  grid  points  but  which  propagate  over  long 
distances . 

The  grid  point  A*  values  are  needed  in  reevaluation  which  follows. 
In  the  first  cycle  of  the  program.  A*  need  be  calculated  only  at  data  points, 
since  the  purpose  of  the  first  cycle  is  to  evaluate  the  reports . 


3.5  Reevaluation  and  Lateral  Rejection 

Knowledge  of  the  blended  pressure  and  resultant  grid-point 
reliability  enables  the  program  to  reevaluate,  individually,  any  information 
that  entered  the  analysis.  By  removing  an  individual  datum  value  and  its 
reliability  from  a  grid  point  and  comparing  it  with  what  remains,  or  the 
"background",  the  true  reliability  of  the  data  may  be  quantitatively  assessed 
(Section  2.2). 

Removal  of  a  report  from  the  blended  fields  is  the  exact  reversal  of 
assembly: 

Ab  =  A*  -  Aa  ,  (69) 

PB  -  (a*  P*  -  Aa  Pn)/AB  ,  (70) 
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where  p  is  the  report  that  was  assembled  with  weight  A  and  p„  is  the 
n  a 

background  pressure,  with  weight  A  ,  that  remains  after  removal  of 

B 

the  report.  The  actual  difference  of  these  estimates  is  p  -  p„  while  the 

-1  -1  1/2  n  B 

expected  or  standard  difference  is  a  =  (A  +  A„  )  ,  where  A  is  the 

n,B  n  B  n 

true  weight  of  the  total  report.  (Expected  difference  must  be  based  on 

true  report  weight;  the  assembled  weight,  A  ,  is  reduced  from  A  when 

a  n 

different  values  are  reported  from  the  same  station,  and  when  a  report  is 
off-time.)  By  squaring  the  actual  and  expected  differences,  the  ratio  may 
be  expressed  in  terms  of  variance: 


np 


ab  a 


B 


TT  (p„-pB)‘ 

n  ' 


(71) 


Some  allowance  may  be  made  for  disparity  in  position  between  p  and  p  , 

n  B 

accounting  for  the  fact  that  reports  occurring  over  a  finite  area  must  be 

expressed  in  a  single  grid-point  value.  This  allowable  disparity  should  be 

2 

a  function  only  of  the  gradient  in  the  area  of  a  grid  point.  If  the  A 
equation  is  written  in  terms  of  blended  values  only,  this  allowance  may 

^  2  •fc 

be  made  by  reducing  the  difference  (p  -  p  )  for  a  portion  of  the  p 

*  n  2 

gradient.  Substituting  p  and  A  in  the  A  equation: 


np 


(A*-A>„ 


A  -  A  +  A 
a  n 


Pn  - 


A*  *  A  1  2 

A  P  "A*Pn 

a  n 


A  *  A 

A  -  A 


a  J 


,  or 


np 


A  A 
n 


*2 


(pn-p*) 


(a*  -  A  +  A  Ha*  -  A 


0(a*-aJ 


(72) 


The  mean  square  gradient  per  unit  grid  length  around  grid  point 
(X,w)  is 

~2  /  *  2  *2  *2  *  2\ /2 

gl,m  =  y'l.m  +Vl-l,m  ) 


(73) 
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A  reasonable  position  disparity  to  allow  is  the  r.m.s.  distance  of  a  large 

2  l/2 

population  of  reports  from  their  respective  grid  points,  or  (d  /6)  , 

where  d  is  the  grid  spacing.  Then  g 2  should  be  weighted  by  the  mean 

ztm 

square  distance  in  unit  grid  spacing,  or  1/6. 

2 

The  final  equation  for  A  becomes 

np 


(74) 


The  value  of  A  measures  how  accurate  a  pressure  report  actually  is, 

as  compared  to  its  expected  accuracy  given  by  the  reliability  assigned  to  it. 
.  2 

The  median  A  for  all  reports  in  the  analysis  should  therefore  be 
np 

approximately  1.0,  indicating  a  balance  between  reports  with  more  error 
and  reports  with  less  error  than  expected.  This  offers  a  method  for 
checking  the  assignment  of  class  reliabilities  to  reports,  which  otherwise 
would  be  arbitrary. 


The  value  of  A  for  a  report  can  now  be  used  for  reevaluation. 

2 

If  A  ^  1  the  report  is  within  its  expected  error  and  no  change  is  made 
np  2  2  2 

in  its  reliability.  If  1  <  A  <  A  ,  where  A  is  the  maximum 
2  np  maxp  maxp 

allowed  A  for  a  pressure  report,  the  report  is  reevaluated  by  the  formula 


2  A 
_ n 


1  + 


(75) 


2  2 

where  R  indicates  reevaluation.  Any  report  of  A  >  A  is  withheld 

np  maxp 

from  the  final  analysis  cycle.  Since  this  rejection  is  based  on  current 
information  in  the  vicinity  of  the  report,  it  may  be  termed  a  "lateral" 
reject  as  distinguished  from  the  initial  "gross  error"  reject. 
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Since  blended  difference  reliabilities  (B*  ,  C*)  are  not  produced, 

the  reevaluation  of  wind  reports  must  be  done  in  a  different  manner.  The 

2  2 

expected  variance  of  a  wind  report  is  approximately  l/B  +  a,  where  o 

n  b  b 

is  the  variance  added  at  the  grid  point  for  the  balance  approximation. 

(This  is  true  only  for  a  single  wind  report  at  a  grid  point  since  the  balance 
variance  is  added  only  once,  independent  of  the  number  of  reports.)  The 
variance  ratio  for  a  wind  report  is 


(„n-,*)2  +  (V„-V*)2 

1/Bn  +  °b 


2  2 

If  1  <  X  <  X  the  report  is  reevaluated: 

nw  maxw 


BnR  °nR 


(1/Bn  +  ffb)  (1+Xnw)/2-*b 


The  subtraction  of  is  necessary  because  it  will  be  added  back  in  the 

final  assembly.  As  in  the  case  of  a  pressure  report,  no  change  in 

2  '  2  2 
reliability  is  made  if  X  SI  and  the  report  is  rejected  if  X  >  X 

nw  nw  maxw 

A  vector  sum  and  difference  ratio  is  also  computed: 


0Jn-O2  +  (t'n-t'*)2 

c+(lln  +  l/:)2  +  (<'n+v*) 


Using  a  value  for  c  of  0.4,  the  report  is  rejected  if  r  >  0.5.  This  test 
is  most  severe  in  light  wind  speed  cases  and  causes  rejection  of  many 
light  or  calm  nighttime  wind  reports  which  occur  as  a  result  of  high 
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boundary  layer  stability.  High  elevation  and  valley  reports,  which  are 
also  nonrepresentative  of  the  sea-level  pressure  gradient,  are  rejected  by 
this  test  as  well. 

The  relative  merit  of  the  separate  first-guess  fields  and  their 
derivatives  can  also  be  evaluated.  In  the  case  of  the  first-guess  pressure, 
the  ratio 


(p,-p*)2 

F  ■  (p6-p*)2+  (px-p*)2 

is  computed  where  p  is  the  prognostic  field  and  p  is  the  extrapolated 

6  X 

field.  The  separate  first  guesses  can  now  be  reweighted  to  form  the  new 
first  guess: 


poR  -  ppx  +  (1  -  0  pe  •  (80) 

The  first-guess  reliability,  Aq  ,  is  not  changed.  The  same  procedure  is 
applied  to  the  first-guess  difference  and  Laplacian  fields.  The  new  first- 
guess  fields  and  the  reevaluated  data  are  saved  for  the  reanalysis  cycle. 


3 . 6  Final  Analysis 

The  purpose  of  the  first  analysis  cycle  is  only  to  get  an  accurate 
evaluation  of  all  data  and  first-guess  fields  for  use  in  the  final  analysis. 

The  reanalysis  begins  by  returning  to  the  assembly  stage  of  the 
program.  The  new  assembly  is  initialized  with  the  reweighted  first-guess 
fields  and  the  reports  are  assembled  with  their  reevaluated  weights.  Gross 
error  and  lateral  rejects  are  omitted. 
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The  blending  for  p*  follows  exactly  the  same  procedure  described 
in  Section  3.3  except  that,  to  aid  convergence,  the  iteration  is  initialized 
with  the  p*  field  computed  in  the  first  cycle.  The  A*  values  are  not 
needed  by  the  program  in  the  final  analysis  but  may  be  computed  for  product 
use.  The  program  skips  the  reevaluation  stage  and  proceeds  to  the  output 
section.  The  final  analysis  is  stored  on  the  disc  in  125x125,  89x89  and 
63x63  grid-point  form.  The  smaller  arrays  are  used  as  input  to  other 
FNWC  programs.  The  full  resolution  analysis  is  saved  for  future  first- 
guess  use  by  FIB/SLP. 


4 .  Sample  Analysis  and  Evaluation 

A  section  of  a  FIB/SLP  analysis  is  shown  in  Fig.  11.  Further 
evaluation  of  FIB/SLP  should  consider  the  character  shown  in  the  analyzed 
fields,  discrimination  in  rejecting  reports,  and  use  of  first-guess 
information.  A  complete  evaluation  of  FIB/SLP  will  be  possible  only 
after  a  suitable  period  of  operational  testing  and  use. 


-65- 


Appendix 


FIB/SLP  Adjustable  Constants 


Symbol 

Value 

Purpose 

f 

o 

16  mb2 

computation  of  D 

o 

(constant) 

£i 

1.2 

computation  of  Dq 

(L2  term) 

f2 

4  mb 

pressure  gross  error  check 

(constant) 

f3 

2.8 

pressure  gross  error  check 

(  |  g  |  term) 

f4 

20  mb 

pressure  gross  error  check 

(1  -  cos  <pterm) 

f5 

0.5  mb2 

computation  of  A^ 

(constant — minimum 
report  variance) 

f6 

0,003 

computation  of  A^ 

(g2  term) 

f7 

,  -6  ,2  -2 
1x10  mb  m 

computation  of  A^ 

(E2  term) 

f8 

2  -2 

0 . 1  mb  hr 

computation  of  A^ 

/  2 

(H  — constant  term) 

f9 

0 . 02  hr  2 

computation  of  A^ 

(H2g2  term) 

£10 

1.0  mb2 

computation  of  B  ,  C 
n  n 

(constant — min .  variance) 

x2 

maxp 

15 

lateral  reject  limit,  pressure 

x2 

maxw 

8 

lateral  reject  limit,  wind 
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